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Abstract 

This paper addresses some theoretical modelling and control issues for a mixing chamber used in 
rocket engine testing at NASA Stennis Space Center. The mixer is responsible for combining high 
pressure LH2 and GH2 to produce a hydrogen flow that meets certain thermodynamic properties 
before it is fed into a test article. The desired properties are maintained by precise control of 
the LH2 and GH2 flows. The mixer is modelled as a general multi-flow lumped volume for single 
constituent fluids using density and internal energy as states. The set of nonlinear differential 
equations is modelled in the SIMULINK environment including a table look-up feature of the fluid 
thermodynamic properties. A small-signal (linear) model is developed based on the nonlinear model 
and simulated as well. Pulse disturbances are introduced to the valve positions and the quality of the 
linear model is ascertained by comparing its behavior against the nonlinear model simulations. Valve 
control strategies that simulate an operator-in-the-loop scenario are then explored demonstrating 
the need for automatic feedback control. Finally, classical optimal single-output and multi-output 
Proportional/Integral controllers are designed based on the linear model and applied to the nonlinear 
model with excellent results to track simultaneous, constant setpoint changes in desired exit flow, 
exit temperature, and mixer pressure, as well as to reject unmeasurable but bounded additive step 
perturbations in the valve positions. 


1 Background 

The Engineering Division at NASA John C. Stennis Space Center (SSC) continues its efforts to 
assemble a software simulation package that captures the static and dynamic characteristics of modem 
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and future thermodynamic systems. The package is foreseen to fulfill the need for an accurate and 
verifiable thermodynamic system simulation environment [1, 2, 3, 4], 

In this work, we build mainly upon [4] and focus on the point-to-point tracking control problem 
for an LH 2 and GH 2 mixer subsystem and associated flow controllers. The relevant components of 
the mixer subsystem are shown in the schematic of Figure 1. The LH 2 valve controls the flow of high 
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Figure 1: Diagram of HPH2 Mixer and Control Valves. 

pressure liquid H 2 from a pressurized tank; the GH 2 valve controls the flow of gaseous H 2 from a set of 
high pressure bottles; and the output valve controls the flow of H 2 into the test article. The primary 
objective of the mixer and control valves is to regulate the thermodynamic properties of the out flow of 
H 2 for optimal performance in spite of the steady depletion of the source bottles, measurement errors, 
modelling imperfections, and other uncertainties and perturbations. 

In [4] , we developed an analytical model of the mixer and control valves suitable for dynamic analysis 
and control design in the vicinity of a thermodynamic equilibrium point. An optimal linear-quadratic 
regulator was designed and simulated using Matlab to illustrate the effectiveness of such a controller in 
compensating against perturbations that result in small deviations around the desired equilibrium. The 
simulations also indicate how one could study which valve or combination of valves is more effective in 
achieving a desired performance (control authority or relegation). 

In this article, the nonlinear mixer equations are simulated in the SIMULINK environment, including 
a NIST-12 property data look-up table incorporated in the model and continuously accessed via Matlab 
functions [5]. The small-signal model developed in [4] requires the numerical computation of several 
partial derivative terms which in turn depend on the equilibrium point. Since there are no analytical 
expressions that would allow the direct calculation of some of these partial derivatives, a Matlab routine, 
in conjunction with the property look-up table, computes these terms as needed. The model is then 
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used to evaluate open-loop control strategies simulating an operator-in-the-loop scenario and testing that 
reasonable results are obtained based on past operator experience. In the future, we hope to validate 
the model with real test data as it becomes available. The simulations section describe closed-loop 
control strategies in more detail based on 1) an optimal multi-input/single-output (MISO) PI controller 
to track the exit flow only; and 2) an optimal multi-input/multi-output (MIMO) PI controller to track 
the exit flow, the mixer temperature, and the exit temperature simultaneously and in the presence of 
step perturbations in the valve positions. 

2 Mixer and Valve Models 

The mixer has been modelled as a general multi-flow lumped volume for single constituent fluids 
using density and internal energy as states [1, 4], The model was independently derived in [5]. In this 
article, we refer to the mixer’s internal energy and density as the “state”, not to be confused with the 
definition of state of a system as described by its thermodynamic properties. Letting Z\{t), Zi{t) denote 
the mixer’s specific internal energy and density, respectively, the two dynamic equations for a constant 
volume V are given by 

Zl — y~- j E net — W ne tZi + <0] (1) 

i 2 = yWnet (2) 

where E net is the net (input minus output) energy into the mixer, W net is the net (input minus output) 
flow, and Q is the total heat transfer rate. The net flow is simply 

W^t — W g + Wi - W e 

where the subscripts stand for “gas”, “liquid”, and “exit”, respectively. In this study, only positive 
flows are considered, that is, the gas and liquid flows are into the mixer and the output flow is out of 
the mixer. In essence, normal operation requires that each valve always experience a positive pressure 
difference across it thus disallowing flow reversal. The net input energy term is given by 

Enet hgXOg T hlXU[ hyVJg 

where h g , hi, h v denote the enthalpy of the GHz stream, LH 2 stream, and volume block (mixer), 
respectively. 

The liquid and exit valves are modelled by an incompressible fluid flow expression that compute the 
flow Wf given the inlet Pi and outlet P 0 pressures, the density of the inlet stream pi, and the valve flow 
coefficient C v f, as follows: 

Wf — y Pi Pw ( Pq ) G v f = otfC v f (3) 


3 


where p w is the density of water (62.4 lbm/ft 3 ), and we have assumed that P* > P a . The gas valve on 
the other hand is more appropriately modelled by the compressible fluid flow expression 


Wf 


2.857 x 1C T 2 ^ y/P?-P*C vf , P„<Pi< 2 P 0 - 
2.423 x 10 - 2 VTiPiC vf , Pi > 2 P 0 . 


(4) 


where the mass flow is in Ibm/s , the pressure in psi , the density in lbm/ft 3 , and the temperature T in 
°R. By design and under normal operation, the properties of the inlet hydrogen gas and liquid are such 
that Pi, P g > P 0 = P v and also P v > P e ; thus, the possibility of reverse flow is not considered in this 
study but can certainly be added if so required. 


2.1 Small-Signal Model 

Substituting the expressions for the valve flows and the net internal energy and net flow into equa- 
tions (l)-(2), and treating the valve coefficients as inputs to be regulated, the dynamic model of the 
mixer is a 2-state, 3-input system of nonlinear differential equations of the form 


— F 1 (Z 1> Z%, Cvgi CvU C've) 
Z2 — F 2 (zi, Z2 , C V g , C v i, C ve ) 


(5) 

( 6 ) 


where Pi(-) and F 2 (-) are nonlinear functions of the state and valve coefficients. 

For the remainder of this article, we will denote constant or equilibrium values of any variable by 
an upper bar (•). Given constant values of valve flow coefficients C v = [ C vg C v i C ve ] T (superscript T 
denotes transposition) and constant fluid properties, the state of the model 


z{t) = 


Z\ : Internal Energy 

z 2 : Density 


reaches a constant equilibrium point z. Next, consider perturbing such an equilibrium by small signals 
x(t) and u(t) so that 


z(t) — z + x(t) 


and 


C v (t) = C v + u(t ) 


where u(t ) denotes a small, valve-coefficient correction/regulation signal. Then, a standard linearization 
of equations' (5)- (6) results in the small signal model 


x 


Ax + Bu 


(7) 


where x{t) is the small perturbation state vector, u(t) is the small perturbation control signal, and the 
two-by-two matrix A and two-by-three matrix B are given by 


,4 = 


dFi 

dFi 


r dFi 

dFj 

dF j 

dz\ 

dZ2 

; B = 

dC V g 

dC vl 

dC ve 

dF 2 

dF 2 

dFz 


dF : g 

dz\ 

dZ2 


dCvg 

dC v i 

dC V e J 
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where the partial derivatives are evaluated at the equilibrium state z and constant valve flow coefficient 
vector C v . 

Output equations of the form 


y = Cx + Du 


where matrices C and D are appropriately dimensioned are easily appended to the model (7) to account 
for the measurement of certain variables such as temperature, pressure, or flow. An output of interest 
is the exit flow w e . Using equation (3) and linearizing around a chosen equilibrium point, we obtain the 
linear approximation 


W elin 


da e da e 


dzi dz2 


x + [a e ] 


eq 


u v 


— C'vfloW'E “ 1 “ flow'U'V 


eq 


A second output of interest is the mixer pressure P v . Taking this pressure to be a function of the mixer 
internal states, that is, 


Pv = fp(z 1 (t),z 2 (t)) 

then, the linearization of this function around the chosen equilibrium gives 
Pvlin ~ 


OR, 8P„ 


dzi dz 2 


X = CpvX 


eq 


A third output of interest is the exit temperature T e . Taking this temperature to be a function of the 
mixer internal states, that is, 


% = Mzi(t),z 2 (t)) 


then, the linearization of this function around the chosen equilibrium gives 


Telin — 


dzi 



X = CteX 


The indicated partial derivative terms are calculated numerically using Matlab functions within the 
SIMULINK model. 

It is not reasonable to assume that the mixer states be measurable. However, the state postulate 1 
implies that the state z(t) — z + x is directly available for control design if at least two independent 
thermodynamic properties were measured. Therefore, 


• assume that real-time measurements of temperature T e and pressure P e are available at the exit 
point. Then, using the thermodynamic table look-up, one finds the exit flow enthalpy h e . 

1 For a simple compressible substance, any thermodynamic property is at most a function of two other independent 
thermodynamic properties. 
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• Assuming a lossless valve, then enthalpy across the valve is preserved so that h v = h e , where h v is 
the mixer (volume) enthalpy. 

• Using h v and T v (measured) and the table look-up, one finds the mixer pressure P v . 

• From h v and P v and the table look-up, one finds the mixer states (energy, density). Or, using 
the measured exit flow w e , the known exit pressure P e required by the specific test article, the 
computed P v , and the exit valve coefficient, one solves Equation (3) for the density p v . 

In summary then, it is reasonable to assume that 


y = x — z — z 


(8) 


is available or that it can be deduced from other real-time measurements. This may be called a “static 
state observer” . An alternative is to design a dynamic observer to provide real-time, optimal estimates 
of the states. This is beyond the scope of this article. 

The calculation of the equilibrium points and the elements of the matrices A, B, C, and D are not 
included due to space limitations. The relevant results are compiled in Appendix A. 

The evaluation of ^ is accomplished numerically via a Matlab routine. The evaluation of ^ on 
the other hand, is done by differentiating the expression 

P 

h = E + C un it s 

P 

where C un i ts = 1 in the SI units or C un i ts = 0.1852 when h, E are in BTU/lbm , P in psi, and p in 
Ibm/ ft 3 . This results in 


Sh'o , ( Cunits S P'o 


dh v 
dz o 


= c, 


units 


22 dzi 


'1 dP v 

Pv 

Z 2 dz 2 

4 J 


It was verified that the relations between pressure P and internal energy or density are essentially linear 
in the vicinity of an equilibrium point. The Matlab functions extract the slope of a line constructed in 
the vicinity of a desired equilibrium or linearization point. 


2.2 A Numerical Example 

It is to be noted that the thermodynamic tables used by the Matlab routines were generated with 
the program MIPROPS and therefore contain an inherent error because of the resolution used. For 
example, the step size in the pressure intervals is 100 psi. However, linear interpolation is used to 
calculate values not found in the tables which reduces the error. Table 1 lists the data for a typical 
equilibrium point EQi corresponding to a given set of outflow requirements. 
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The matrices A and B are evaluated at the equilibrium EQi listed in Table 1 to be 

332.74 -42.01 

0.6612 -0.4919 

[1.23] u ve 

In the absence of real data, we shall refer to the nonlinear SIMULINK simulation as the real system. 

Several simulations were performed to evaluate the linear model: 

• Disturbance Effect - No Control. We ran the model under a 20 percent pulse disturbance in the gas 
and liquid valves that take the state away from the initial equilibrium point. When the disturbance 
disappears, the mixer returns to the equilibrium. It was observed that the linear model captures 
the important dynamic characteristics of the system very accurately. 

• Tracking and Disturbance Rejection - Open Loop Control. The behavior of the mixer model was 
evaluated as the valves are steered in an open-loop fashion from the initial to the final pre-calculated 
positions. These were generated to simulate the manner in which an operator might open/close the 
valves in response to a visual monitoring of pressure, temperature, or flow sensors. The simulations 
included a 20 percent pulse perturbations on the liquid and gas valve positions. 

3 Point-to-Point Tracking Control Design 

The mixer- valve system is a multi-input, multi-output system with an associated transfer function for 
each input-output pair. It is in general difficult to design JPID-style controllers for each loop because 
of the interactions between channels. The approach chosen in this work is based on the state-space 
model (7) rather than on the individual transfer functions. 


A = 


-17.52 -3393.08 

-0.1597 -43.45 


B = 


25.11 

0.651 


For an output y = w e , the linear approximation is 

'Melin C v fi ow X -f- D flow'U've = [0.2386 66.84] X + 


3.1 Single Output Tracking 

Consider the problem of exit flow tracking and define a tracking error 

e(tj — w e (t) (9) 


where the reference or desired constant exit flow is w er . We seek to design a controller that steers the 
error signal to zero, that is, e(t) — * 0 as t — » oo. To that end, define an augmented state vector 


X = 


x 

e 


( 10 ) 
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that satisfies the differential equation 


A 

0 ’ 


B 



X + 

ii = AX + Bit 

-C flow 

0 

0 0 — Dfiow 



Given the system (11), a performance criterion of the form 

J = J (x T (t)QX(t) + u T (t)Ru(t)^j dt 

is minimized where the matrix Q = Q T > 0 and the matrix R = R T > 0. The needed controllability 
conditions for the pair (A, B) are by now standard and may be found in [6]. The resulting optimal 
controller is given by the linear full-state feedback law 

u(t) = — KX (t) (12) 

where the matrix K is constant and is found a-priori by solving the steady-state Riccati Equation 


A T P + PA- PBR^B? P + Q = 0 =* K = R~ 1 B T P 


The closed-loop system 

x = [A-bk]x 

is guaranteed to be asymptotically stable, that is, X(t) — *■ 0, t — » oo. Note that driving X to zero is 
equivalent to driving the error e to zero as required, and the mixer state a: to a constant. Finally, integrat- 
ing both sides of (12) leads to the optimal multi-input/single-output proportional-integral (MISO-PI) 
controller expression 

u{t) = - K P x(t) -Kj [w er - We(r)\ dr (13) 

J o 

For example, select the following weighting matrices: 


■ i 

0 

0 


' 0.1 

0 

0 

0 

1 

0 

; R = 

0 

1 

0 

. 0 

0 

20 _ 


0 

0 

0.1 _ 


so that it is equally important to steer x\ and to zero, but it is 20 times more important to drive the 
exit flow error to zero. Similarly, the liquid and exit valves have 10 times more authority than the gas 
valve in achieving the stated goal. The optimal gain matrix is found to be 


K - 


0.4679 22.5556 -12.1274 

0.8992 -9.3271 -0.2184 

-0.8653 13.2205 -7.2421 
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from which 



0.4679 

22.5556 1 


' -12.1274 ' 

K P = 

0.8992 

-9.3271 

; Kj = 

-0.2184 


-0.8653 

13.2205 _ 


-7.2421 


The final equilibrium point was selected to be that corresponding to a 20 percent increase in exit 
flow, from w e = 37 Ibm/s to w e = 44.4 Ibm/s. The resulting valve, mixer states, flows, and linear 
model performance was evaluated by simulations. There are (small) step changes that appear in the 
temperature and pressure plots caused by the resolution of the thermodynamic table and not by the 
controller. The PI controller successfully takes the exit flow to its final required value. 

Simulation results under the same PI controller and with the added perturbations on the liquid 
and gas valves show that the controller recovers from these perturbations. This behavior could not be 
achieved under open-loop schemes. However, the mixer states, and the valves do not go to their unique, 
pre-calculated final values as would be given by a steady-state equilibrium calculation. This is because 
the controller is designed to track the exit flow step reference only. This is corrected in the next section. 


3.2 Multi Output Tracking 

Next, we design an optimal multi-input /multi-output PI (MIMO-PI) controller that simultaneously 
regulates the exit flow, the mixer pressure, and the exit temperature to arbitrary constant set points 
w er , P vr , T er , respectively. The new augmented vector is 


X = 


x 

e 


VJ er 

p 
1 vr 

T 

± f> r 


x 

- W e (t) 
~Pv(t) 

- T e (t) J 


that satisfies the differential equation 


X = 


A 

0 

0 

0 ' 



B 


Cflow 

0 

0 

0 

X + 

0 

0 

T) flow 

CjpV 

0 

0 

0 

0 

0 

0 

-c te 

0 

0 

0 _ 


_ 0 

0 

0 


ii — AX -j- Bii 


(14) 


Minimizing the quadratic cost introduced in the previous section leads to the linear MI MO — PI 
controller 


u(t) = - K P x(t) - Ki f* 

Jo 


w er - w e (r) 
Pvr - Pv(r) 
Ter ~ T e (r ) 


dr 


(15) 
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As a simulation example, select the following weighting matrices: 

'1 0 0 0 O' 

0 1 0 0 0 1 0 0 ' 

Q = 00 20 00 ; R = 010 

00010 [001. 

0 0 0 0 200 _ 

so that it is equally important to steer xi, X 2 , and the mixer pressure error to zero, it is 20 times more 
important to drive the exit flow error to zero, and it is 200 times more important to drive the exit 
temperature error to zero. Several simulations indicated that the temperature tracking is the slowest 
of the three and hence the more stringent requirement. This time, all Valves have equal authority in 
achieving the stated goals. The resulting optimal gain matrix is found to be 

' -0.0719 46.93 -1.8728 -0.8897 2.5733 ' 

K = 0.9454 -9.8472 -0.384 -0.1607 -13.91 

. -0.0439 -11.455 -4.0429 0.4274 0.1296 . 

Figure 2 is the SIMULINK diagram used in the simulations. The details of the mixer are not shown. 
The set point requirements are a 20 percent increase in exit flow from 37 to 44.4 Ibm/s ; a 100 psi 
increase in mixer pressure from 6000 psi ; and a 2 °F increase in exit temperature from —240 °F. These 
are introduced as step commands at time t = 0.3 sec. The valve disturbances are set on and a rate 
limiter with a slew rate at 50 is put in series with each valve to include another level of realism in 
the simulations. The simulation plots are in Figure 3 through 6. The effect of the rate limiter is an 
increase in overshoots and reduced overall performance quality while precluding unrealistically fast valve 
commands rises especially at the start of the command. The outputs are controlled as required. It is 
expected that real valves will not only exhibit rate limiting but other linear and nonlinear characteristics. 

4 Conclusions 

An analytical small-signal model of the high-pressure H 2 mixer and control valves has been derived. 
The model is suitable for a variety of analysis and control designs valid in the vicinity of an operating 
or equilibrium point. A classical approach to control design would involve the consideration of a total 
of nine single-input, single-output transfer functions. Unless the control loops are weakly coupled, 
such a design is difficult to complete because of the interaction among the control channels. Instead, 
a state-variable control design has been pursued that results in an optimal multi-input/multi-output 
proportional/integral law that can be implemented in real-time. It is straightforward to design the 
controller and it is possible to assign different authorities to each control valve. Precise and simultaneous 
tracking of constant setpoint changes in exit temperature, exit flow, and mixer pressure were achieved 
in the presence of constant and unmeasurable additive valve position perturbations. Current efforts 
include the validation of the model with test data. 
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l 



P 

(psi) 

T 

(F) 

P 

(Ibm/ ft 3 ) 

h 

(BTU/lbr n) 

C V open 

%OPEN 

w 

Ibm/s 

21 

Z2 

gh 2 in 

13500 

90 

2.91 

2113.6 

230 

0.96 

- 

- 

~ 

gh 2 out 

6804 

130.4 

1.702 

2113.6 

- 

- 

3.72 

- 

” 

lh 2 in 

8500 

-340 

5.042 

329 

115 

17.74 

- 

- 

- 

lh 2 out 

6804 

-326 

4.519 

329 

- 

- 

33.21 

- 

- 

Mixer 

6804 

-273.32 

3.84 

508 

- 

- 

| 

181.39 

3.84 

Outflow 

5533 

-266.21 

3.404 

508 

270 

11.12 

36.93 

- 

- 


Table 1: Initial Equilibrium Data. Underlined Values are Computed. (IN=Into Valve; OUT— Out of 

Valve), zi in BTU /lbm\ z 2 in Ibm/ ft 3 . 
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Appendix A 

Given a set of required exit mass flow, pressure, and temperature w e , P e , T e , then in steady-state [4] 

(16) 


C V g — 

hy hi 1 _ 
=^—w e 

hg hi CXg 

_ hy hi 

^ Wg ~ hg-h 

II 

hg h v 1 _ 
= w e 

hg - hi ai 

_ hq hy 

Wl ~ hg-h 


W P 


W P 


(17) 


To determine the linear perturbation model (7), it is necessary to find the analytic expressions of the 
indicated partial derivatives. For simplicity we include only those terms that do not evaluate to zero at 
the equilibrium point. Omitting the details, these are found to be: 


dF, 


dz\ 

dFi 

dz 2 

dFj 
dzi . 

~ m 

dz 2 

' m 

_dC V g J 

' dF\ 
dC ve 

dF 2 

dC vl 


eq 


eq 


eq 


eq 


Fz 2 

1 

Vh 


' da q „ dai . . da e dh v 

(hg — Z l) Cvg-Q^- + (hi — Zi) C v l (h v Z\ ) C ve W e 


J e 9 


. _ da„ ,, dai ,, , „ da e dh v 

(hg - Zi) Cvg-g^- + (hi - Zi) - (hv - 2l) Cve g^ w e 


' dz 2 


eq 


1_ 

V 

l_ 

V 


r da 9 

V9 d Zl 


dai 
‘ dz\ 


da, 


n l_l _ n 

u vl Cx 


da g dou dae 

"T* {-'vl o {^ve Q 

y UZo OZ 2 OZ <2 


eq 


eq 


= (*»-*) • 


eq 


eq 


dFr 


dC vt 


eq 


(T 

Vz 2 


(hi - zi) 


eq 


a e ( 7 

- 'N 


' dF 2 " 

1 _ 

Vz 2 ( K 

-zi) 

) 

dC vg \ eq 

- v a 9 

1 

' dF 2 ~ 


1 _ 


v &l ; 

dC ve _ 

eq 

~v a * 



where 


da n 


dzi 

dai 

dzi 

dae 

dz\ 


= < 


eq 


eq 


eq 


.2.857 x 10-^ ^Ts 6 ^ 


0, 


P v < Pg < 2 P v ; 
Pg > 2 P v . 


-1 

1 

2 


pip w dP v 


ai dzi 

pvpw dP v 

a e dz\ 


eq 


eq 


12 


da g 

dz 2 


eq 


-2.857 x 10-^ 

0, 


Py dPy 

y/PS=F2»* 


P v <P g < 2 P v \ 
Pg > 2 Pv. 


da t 

-1 

f* Pipw dP v 


da e ' 

1 

Pw 1 

dz 2 _ 

2 

eq 

8 

n 

to 

1 

5 

eq 

dZ2 m 

2 

e? 

Va e l 


Pe + 



eq 
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Figure 2: SIMULINK Diagram of the Mixer with MIMO PI Control Loop. 
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Figure 3: MIMO PI Control vs Time (sec) (with 
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Figure 4: Mixer Behavior vs Time (sec). MIMO PI Control (with valve perturbations) 
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Figure 5: Input and Output Flows with MIMO PI Control vs Time (sec) (with valve perturbations) 
Flow Balance is w g + wi — w e . 


Exit Flow (NL Model) - Ibm/s 



Exit Temperature (NL Model) - deg F 



Figure 6: Exit Flow and Temperature Error vs Time (sec). (MIMO PI Control with valve perturbations. 
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